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The stability of nonstationary states of homogeneous spin-1 Bose-Einstein condensates is studied 
by performing Bogoliubov analysis in a frame of reference where the state is stationary. In particular, 
the effect of an external magnetic field is examined. It is found that a nonzero magnetic field 
introduces instability in a 23 Na condensate. The wavelengths of this instability can be controlled by 
tuning the strength of the magnetic field. In a 87 Rb condensate this instability is present already at 
zero magnetic field. Furthermore, an analytical bound for the size of a stable condensate is found, 
and a condition for the validity of the single-mode approximation is presented. Realization of the 
system in a toroidal trap is discussed and the full time development is simulated. 
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I. INTRODUCTION 

The excitations and stability of spinor Bose-Einstein 
condensates (BECs) have been the subject of intense 
study in recent years. The topic was first explored in two 
seminal theoretical papers in 1998 0,3- These papers 
discussed the stability of F = 1 BECs against small ex- 
ternal perturbations, such as fluctuations in the trapping 
potential or small magnetic field gradients. In an unsta- 
ble spinor condensate small perturbations may lead to 
domain formation, where the populations of spin compo- 
nents become position dependent. The theoretical stud- 
ies typically examine the stability of stationary states. 
This is done using a linear stability analysis, where a 
small perturbation is added to the stationary state and 
the time evolution equations are expanded to first order 
in the perturbation 0tM| • These theoretical studies have 
shown that the stationary states of ferromagnetic F = 1 
condensates may be unstable. Antiferromagnetic con- 
densates in stationary states, on the other hand, appear 
to be stable against small perturbations in the absence 
of external fields. 

The properties of spinor condensates have been stud- 
iedby various experimental groups (see, for example, 
|9i-tl7|). Signs of instability have been observed by the 
Chapman group fllj ] in a ferromagnetic F = 1 conden- 
sate and the Sengstock group in an antiferromagnetic 
phase of anF = 2 rubidium condensate [l4[ ■ The group 
of Stamper-Kurn has observed that an F = 1 87 Rb con- 
densate prepared in a stationary paramagnetic state de- 
velops spin textures as it is rapidly quenched across a 
quantum phase transition flM [l6j |. They saw a similar 
phenomenon when an unmagnetized rubidium gas was 
cooled to quantum degeneracy [l~7j . Experiments often 
concentrate on spin-mixing dynamics, which can be ini- 
tiated by preparing the condensate in a non-stationary 
state. There are a few theoretical studies where the sta- 
bility analysis has been extended to nonstationary states: 



The effects of a non-zero magnetic field on the stabil- 
ity of states with time-independent spin populations but 
oscillating relative phases have been examined in Rcfs. 
[l8l - |20| . Although the spin populations remain constant 
throughout the time evolution, these states are nonsta- 
tionary because the relative phases of the spin compo- 
nents vary in time. Another example is given by Ref. 
[2l| . where the stability of states that show oscillations 
both in spin populations and relative phases was studied, 
under the assumption that the magnetic field vanishes. 

In the present paper, we generalize these findings to ar- 
bitrary states with and without magnetic fields. We con- 
centrate on the case where the magnetic field is nonzero 
and also consider a situation where the spin populations, 
as well as the phases of the spin components, oscillate 
in time. In the analytical calculations we assume that 
the particle density is homogeneous. We solve analyti- 
cally for the eigenmodes and eigenstates of a certain class 
of nonstationary states in rubidium and sodium conden- 
sates in an arbitrary magnetic field. We show that this 
can be done in a simple way by transforming to a refer- 
ence frame where the states in question are stationary. In 
this way, the time dependence of the matrix determining 
the stability properties can be eliminated, and the prob- 
lem becomes easily tractable. We emphasize that we are 
studying the stability in nonzero magnetic field. The 
zero-field case has been discussed in (21J. Knowing the 
eigenmodes allows us to derive an analytical formula that 
connects the stability of a condensate to its size and the 
strength of the magnetic field. In [l8| it was found that 
a sodium condensate is unstable at a low magnetic field 
provided that the condensate is larger than the spin heal- 
ing length. We show here that in a situation where the 
magnetic-field energy dominates over the spin interaction 
energy, both rubidium and sodium condensates may be 
unstable even when the size of the condensate is smaller 
than the spin healing length. 

This paper is organized as follows. Section [TTJ intro- 
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duces the system and presents the Hamiltonian and time 
evolution equations. Sec. IIIII formulates the theory of 
the Bogoliubov analysis of nonstationary states. In Sec. 
IIVI analytical results concerning the stability of F = 1 
condensates are derived. In Section [V] the stability of a 
state orthogonal to the magnetic field is studied using 
the results of Floquct theory Additionally a sufficient 
condition for the size of a stable condensate is derived 
and a condition for the validity of the single-mode ap- 
proximation is presented. In Sec. I VII the realization of 
instabilities in a condensate confined in a toroidal trap is 
discussed. The time development is simulated using the 
Gross-Pitaevskii equation. Finally Sec. I VIII contains the 
concluding remarks. 



II. THEORY OF A SPIN-1 CONDENSATE 

The order parameter of a spin-1 Bose-Einstein con- 
densate can be written as ip = (ipi,ipo,ip-i) T , where T 
denotes transpose. The normalization is now chosen as 
Sm=-ilV'm| 2 = n i where n is the total particle den- 
sity. We assume that the trap confining the condensate 
is such that all the components of the hypcrfine spin can 
be trapped simultaneously and are degenerate in the ab- 
sence of magnetic field. If the system is exposed to an 
external magnetic field that is parallel to the z axis, the 
energy functional reads 

Ety] = [ dr{^{v)H{v) + Ug [^{v)^{v)] 2 + g 2 (F) 2 } 



p(F z )+q(F 2 )}, 



(1) 



where h = + U(r) - fi , F = (F Xl F y ,F z ) is the 

spin operator of a spin-1 particle and we use the nota- 
tion (X) = ^'(r)X^(r). Here U is the external trap- 
ping potential and the chemical potential, taking care 
of the conservation of the total particle number, is de- 
noted by \i. The strength of the spin-independent inter- 
action is characterized by go = AitY?(ao + 2a2)/3m, while 
#2 = 47rft 2 (a2 — ao)/3m describes the spin-dependent 
scattering. Here ap is the s-wave scattering length for 
two atoms colliding with total angular momentum F. 
For 87 Rb the scattering lengths used in this paper are 
a = 101. 8a,B and a 2 = 100. 4a# [22| with ag being 
the Bohr radius. For 23 Na the corresponding values are 
a = 50.0a_B and a 2 = 55. la^ Q. (Note, however, 
that there are many estimates for the difference a 2 — ao 
in the literature 0, [H, [H, The magnetic field in- 

troduces two terms, one of which is given by the linear 
Zeeman term p = —gp,^B, where g is the Lande hyper- 
fine g factor, /iB = eh/2m c is the Bohr magneton (m e is 
the electron mass, and e > is the elementary charge), 
and B is the external magnetic field. The other term is 
the quadratic Zeeman term 



E hf 



where Em is the hypcrfine splitting. For 87 Rb and 
23 Na the hyperfine splittings are -Em = 6.835 GHz and 
Em = 1.772 GHz, respectively. In both cases g = —1/2. 
The value of q can be made negative by using a linearly 
polarized microwave field 25j. In this paper we concen- 
trate on non-negative q. 

We characterize the spin of the state ip by the spin 
vector f , defined as 



f(r) = 



V> T (r)F^>(r) 
n(r) 



(3) 



The length of this vector is denoted by /, / = ||f||. In 
addition to the number of particles, the magnetization in 
the z direction, defined as 



M, 



J drn(r)f z (r) 
Jdrn(r) 



(4) 



is also a conserved quantity. The Lagrange multiplier 
related to magnetization can be included into p. We con- 
sider mostly homogeneous systems, for which M z = f z . 
The time evolution is governed by 

ih-^(t) = H[m]m, (5) 

where the Hamiltonian is defined as 

H[ip] = [h + 50 n(r)]I + g 2 (F) ■ F - pF z + qF 2 z . (6) 

For a homogeneous system h — > — /i, and the density n 
becomes position independent. Consequently the energy 
of a homogeneous system reads 

Ebp] = -fi + \ (g n + g 2 nf) - pf z + q/n{F 2 z ). (7) 

In the following analysis the time evolution operator 
Uip of the state ip, ip(t) = U^,(t)ip(0), is used frequently. 
This operator can be formally written as 



(8) 



(2) 



where T is a time-ordering operator. Note that the 
Hamiltonian appearing in the exponent depends on the 
state of the system. In some cases can be solved 
analytically, but in general, numerical calculations are 
necessary. In this paper the numerical calculation is 
done by first solving the time evolution of ip, with the 
help of which we get the time-dependent Hamiltonian. 
The columns of the propagator Uj, can then be obtained 
by calculating the time evolution (under H) of the basis 
states (1,0, 0) T , (0,1, 0) T , and (0,0, 1) T . 



III. STABILITY OF NONSTATIONARY STATES 

We study the stability of nonstationary states using 
Bogoliubov analysis. This is done in a basis where the 
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state we arc interested in is time-independent. We define 
a new (time-dependent) basis {| + l) ncw , |0) now , | - l) ncw } 
in terms of the old basis {| + 1), |0), | - 1)} as |i/) now = 



+1,0,-1. Here EL, is defined as in Eq. 



In the new basis, the energy of an arbitrary state 
given by 



is 



£ ncw [0] =E[U^]+ih{<l>\ 



dt u * 



U^), (9) 



and the time evolution of <p can be obtained from the 
equation 



d(j) _ SE ncw [4>] 
dt ~ Stf ' 



(10) 



Equation (0) can be simplified using the equation 



if). 



U 4 , = -t/" 1 if [</>]£/</>• Using Eqs. ©, © 
and (TTUJ) it is then easy to see that the state </> = tp(0) 
does not evolve in time, confirming that ip(0) is a sta- 
tionary state in the new frame. We study the stability 
of V'(O) by replacing ip(0) — > 4>(0) + 5ip in the time evo- 
lution equations obtained from Eq. (fTU|) and expand the 
resulting equations to first order in Sip. The perturbation 
Sijj = (Sipi, Sipo, 5ij}-i) T is assumed to be of the form 



E 

k 



[«i;k(t) e lk - r -«* k (t)e- 



ik-rl 



-1,0,1. 



Straightforward calculation gives the differential equa- 
tion for the time evolution of the perturbations as 



ih 



dt 



( u l;k \ 




/ Ml;k \ 


"0;k 




"0;k 


W-l;k 


= H B 


W-l;k 


Wl;k 


«l;k 


«0;k 




«0;k 


\«-l;k/ 




Wl;k/ 


(f " 






\Y* - 


-X* ' 





where the 3x3 matrices X and Y are defined as 

x = e k+go \m)(m\ 

+92 \ul{t)hm){ul{t)hm\ 

u—x,y,z 

+ 92 E l^{(*)^(t)><[^{(*)^(*)]*|. 



(11) 



(12) 



2J.2 



Eft 



: 2m 



(13) 

(14) 
(15) 



and tp(t) = (t)i()(0). In the rest of the paper we call the 
operator Hg the Bogoliubov matrix. The magnetic field 
dependence appears in the Bogoliubov matrix through 
the magnetic field dependence of U^. The operator Hb 



is typically time-dependent and the time evolution of the 
perturbations is given by the time-ordered integral 



U B (t) = fe- l/h ti dT " B(T) . 



(16) 



In general, both and Ub have to be calculated nu- 
merically. In some special cases it is possible to express 
Ub analytically in terms of a time-independent Bogoli- 
ubov matrix, and the stability can be determined by cal- 
culating the eigenvalues of this matrix. The system is 
unstable if at least one of the eigenvalues of Hb has a 
nonzero complex part. Another case considered in this 
paper is one where the time evolution of Hb is periodic. 
This makes it possible to use Floquet theory to study the 
stability. We first discuss some special cases that allow 
analytical solution, and then proceed to the case where 
Hb is periodic. 



IV. ANALYTICAL RESULTS 

In this section the stability is studied using mainly an- 
alytical means. First we analyze the stability of a sys- 
tem where the spin and magnetic field are parallel in the 
initial state. In the second case we concentrate on the 
stability in the limit of a large magnetic field. 



A. Parallel spin and magnetic field 

One case where the stability can be studied analyti- 
cally is a system where the spin and magnetic field are 
parallel in the initial state, (F x ) = (F y ) =0. It is easy 
to show that the state has to be of the form 



(i'li 



V(l + /*)/2 N 
_0 



\M=f, (17) 



where the relative phase of the two nonzero spin compo- 
nents can be chosen to be zero due to the fact that the 
energy is invariant under rotations around the z axis. In 
general, the stability properties of two states that can be 
obtained from each other using an element of the sym- 
metry group of the energy are identical [l[. Therefore, 
instead of studying the stability of all possible states, it is 
enough to concentrate on those states that cannot be con- 
nected by an element of the symmetry group. We remark 
that the stability analysis presented in this section is valid 
for all states at zero magnetic field. In this case we can 
make use of the fact that for any spin state ip there exists 
a spin rotation operator R(a, (3, 7) = e~ laF * e~ l P Fy e~ i ~ /Fz 
such that tp = e tT R(a,(3, ^y)ip\\, where (a, f), 7) are the 
Euler angles and r is the global phase. At zero magnetic 
field the initial state can therefore always be assumed to 
be of the form ipu given in Eq. (|17p . 

In Appendix [A] we show that the spin populations of 
ip\\(t) = U^j),, (£)V'|| are time independent regardless of the 
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value of q. Then Eq. (j6)) gives the propagator 

and the matrices appearing in the Bogoliubov Hamilto- 
nian become 

I" =efci + so|^|(0))^|,(0)| 

+ g 2 \F z ^ (0)) (0)| + 5 2 n(I - F 2 ) (19) 

y" = . 9 ol^ii(o))(V||(o)| + 52 |iW||(o)><i^n(o)| 

(20) 



-i2qt/h 



g 2 n^l-f 2 (l - F 2 




i z e* 



i L 



i z <* 



The Bogoliubov matrix Hg is such that the time evolu- 
tion of {wo;k, ^0;k} is decoupled from the time evolution of 
{ui ; kj U-i;k> Wi ; k, «_i ; k}- Moreover, the Bogoliubov ma- 
trix giving the time evolution of {wi ; k, W-i;k, wi ; k, «-i;k} 
is time-independent and the time-dependence of the 
{wo ; k,^0;k} -part can be eliminated by defining a new 
basis «o;k = e 19 *' ft uo ; k, «o,-k = e _t9 *' ft uo ; k) and -Sj ; k = 
%;k,%;k = v.j;k f° r J = il- After this the eigenvalues 
can be easily calculated 



(fiwi,2) 2 = £fe (5o + 32)n + e k + ns/ (g Q - g 2 ) 2 + 4g g 2 / z 2 

(21) 

(huj 3A ) 2 = e k (g + g 2 )n + e k - n^J (g - g 2 ) 2 + %)52/f 

(22)" 



(t^ 5 , e ) 2 = (ff 2 n) 2 (/ 2 - 1) + (e k + g 2 n - qf 



(23) 



For = the eigenvalues (|21 j) -(|23 |) reduce to those given 
in [21(. We assume that g > and \g 2 \ <C 50, which is 
the case both for rubidium and sodium. Now uq^ are 
always real, but uj^a can be complex if g 2 < 0; the un- 
stable states lie inside a triangular region in the (e k , /f ) 
plane; see Figs. nja)-QJc). F° r fixed values of cj 5 , 6 and 
q, equation (|23[) determines an ellipsoid in the (e k ,f z ) 
plane. The unstable states lie in the interior of the ellip- 
soid obtained by setting ui^^ = and are constrained by 
the inequalities e k , f 2 > 0; see Figs. QJa)-(c). For g 2 > 
the region of instability is shifted by 2g 2 n with respect 
to that of the g 2 < system, as can be seen from Fig. 
[TJ We see that ip\\ is unstable in a rubidium condensate 
if \fz\ < 1- The same applies in a sodium condensate 
if > g 2 n. When q < g 2 n, this state is unstable if 
f 2 < —q 2 + 2q. At \f z \ = 1 the system is stabilized by 
the conservation of magnetization. 

Regardless of the sign of g 2 , the fastest-growing unsta- 
ble mode is located at e k = q — g 2 n and corresponds to 
the wavelength 



A 



2irh 



y / 2m(q - g 2 n) 



(24) 



For a sodium condensate in a magnetic field q < g 2 n 
the fastest-growing mode is at e k = 0. In Fig. [2] we 
show the possible wavelengths of unstable perturbations 



FIG. 1. (Color online) The amplitude of the unstable fre- 
quencies ui = Im[w] in the case f || B for (a)-(c) rubidium 
and (d)-(f) sodium. The units of e k and Wi are \g2\n and 
\g2\n/h, respectively. Here in (a) and (d) q — 0, in (b) and 
(e) q — \ga\n, and in (c) and (f) q = 2|(/2|n. The green color 
(left lobes in top row) indicates the unstable modes given by 
Eq. (|22[) [here called magnetization modes; see (|26[) ]. while 
the blue color gives the instability arising from the modes of 
Eq. (|23[) [now called spin modes; see l|27[l]. The region corre- 
sponding to / = in the bottom row agrees with the results 
presented in Fig. 4(b) of Ref. 




(b) 
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FIG. 2. The wavelengths of the unstable perturbations in 
the case B || f for (a) rubidium and (b) sodium. Horizon- 
tal (vertical) lines denote magnetization (spin) modes. We 
have chosen n = 4 x 10 14 cm -3 and / = 0. The latter choice 
gives the largest possible interval of unstable wavelengths. 
The shaded region gives condensate sizes, which correspond 
to stable systems regardless of the initial state; see Sec. [V] 



as a function of the magnetic field. We have chosen n — 
4 x 10 14 cm~ 3 . 

The eigenvectors {x j; k} corresponding to the eigenval- 
ues ([2"Tj) - (|2"3"|) can be calculated analytically and are given 
in Appendix [B] Using the analytical expressions for the 
eigenvectors the corresponding spin states can be calcu- 
lated straightforwardly; see Eqs. (|B1[) and (|B2[) . We 
denote by dtp 1 the state corresponding to eigenvector i. 
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Wc find to lowest order in g 2 / go (see Appendix [Bj 



Sip 



1,2 



Sip 



3,4 





/ 1 / 1 1 f \ 1 c\ \ 

( V(l + .^)/2\ 








W(l-/.)/2/ 




^ v/(l-/,)/2 








W(l + /»)/2. 



5V 



5,6 



£C 5 - 6 (k-r,g,t) ( l) , 



(25) 
(26) 
(27) 



where C^^ +1 contain all position, time, and magnetic 
field dependence. Of these, Sip 1 ' 2 corresponds to a change 
in density, while the magnetization, defined as in Eq. Q , 
and spin direction remain unchanged. We therefore call 
it a density mode. The perturbations (5-0 3,4 , now called 
magnetization modes, affect the density and magnetiza- 
tion but not the spin direction. Finally, Sip 5,6 change 
the density and spin direction but not the magnetization 
and are called spin modes. The density modes are al- 
ways stable, reflecting the fact that the spin-independent 
interaction is now repulsive. For 172 < the magneti- 
zation mode can be unstable, whereas for 172 > it is 
always stable. This can be understood by looking at how 
the energy behaves when the system breaks into regions 
with different spin values. Neglecting constant terms, the 
energy of an arbitrary state can be written as 



E=^g 2 nf 2 + q(l-p ). 
In the initial state ipn the energy reads 



E\ 



I = \92nfj + q. 



(28) 



(29) 



Assume that in a region of length L\ (L 2 ) the expectation 
value of the spin in the z direction is f z \ {f z2 )- The length 
of the spin vector in the (x, y) plane is denoted by f±i and 
f±2- We choose po = as the magnetization modes do 
not populate the zero component. Consequently, = 
/_i_2 = 0, and taking into account the conservation of 
magnetization, we obtain the equations 



1 Ltfa + Lifa , 

E = — Q 2 n — — + Q. 



::2 



Li + L 2 



(30) 
(31) 



Without loss of generality, we choose f z > 0, f z \ > f z , 
and f z2 < f z . With the help of Eqs. ^ and flST} we ob- 
tain E = g 2 nf 2 {xi+x 2 -x 1 x 2 )/2 + q, where Xi = f zi /f z . 
Taking into account that X\>1 and x 2 < 1, we find that 
x\ +x 2 — x\x 2 > 1. Hence, for rubidium E < En and do- 
main formation is energetically allowed. Conversely, for 
sodium E > E^ and region formation is forbidden for en- 
ergetic reasons. Here we have neglected the contribution 
from the kinetic energy. The energy cost caused by the 



kinetic energy allows only structures with long enough 
wavelength compared to the energy gained from the in- 
teraction energy. When / « 1 in the initial state of a 
rubidium condensate, this energy gain is very small and 
allows only structures with a very long wavelength. This 
qualitative result agrees with Fig. QJa)-(c). 

The spin mode (|23|) increases the population of the 
zero component. Hence we assume domains such that 
Po = 1 = fzi = 0) and p Q = f ±2 = 0,f z2 = 1. As 
before, we have also chosen f z > 0. We get 



( \gtnf z - q 



(32) 



For rubidium this is negative regardless of the value of q, 
and domain formation is possible. For sodium the mag- 
netic field has to be nonzero for instability to appear. As 
q increases, the energy difference E — En grows. This 
excess energy is transferred into kinetic energy of the do- 
main structure. For large enough q this kinetic energy 
has a finite minimum value, and consequently, the wave- 
lengths of the unstable perturbations are bounded from 
above 0. This is illustrated by Figs. [TJb) andQJc). 



B. Stability when q ^> \g2\n,tk 

Another case where it is possible to obtain analytical 
results concerning the stability of the system is when q 3> 
Ifl2|ft, £fc- The relevant parameters characterizing the spin 
states can be determined by writing the general spin state 
as 



gen 



f~ VT —iaF z - 

*ne e e 







(33) 



Here /3 gives the angle between the z axis and the spin 
direction. The global phase r is irrelevant and will be set 
to zero. Furthermore, due to the invariance of the energy 
in rotations around the z-axis, we can choose a = 0. The 
important parameters are then /3 and 7. In Appendix IA1 
we derive an approximate propagator for the system in 
the limit q 3> |<?2|w. It is given, up to a time-dependent 
phase, by 



-it[(g2ncos/3f- P )F z + (g 2 n(2p -l)+q)F^]/h 



(34) 



where po is the initial population of the |0) component. 
When analyzing the stability as a function of /? and 7, 
it is important to note that fixing the direction of the 
spin does not fix the populations: The spin direction is 
determined by /3, while 7 controls the populations of the 
spin components. In more detail, 



Po 



= i[l-V / W I cos(2 7 )]sin 2 /3. 



(35) 



Now (3 and 7 can be chosen to lie in the interval [0, tv/2] 
as the stability properties are identical for states corre- 
sponding to p and 7r — p and similarly for 7. For fixed /3 
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and /, the population po is minimized (maximized) when 
7 = (7 = tt/2). 

Using Eqs. ([I~3" |) . (fTl)) , and ([M]) we obtain a Bogoliubov 
matrix where the time dependence appears via terms of 
the form e ±2l 9*/^ We use the rotating wave approxima- 
tion and set these terms equal to zero. This approxima- 
tion can be assumed to be valid when the quadratic Zee- 
man term is much larger than the other energy scales, 
q 3> efc , I <7a I . The eigenvalues of the resulting time- 
independent matrix can be calculated analytically, but 
they will not be presented here as they have a very com- 
plicated form. The eigenvalues show that a sodium con- 
densate is always stable against long wavelength pertur- 
bations, which is in agreement with the results of the 
previous subsection if q gin. Rubidium condensate 
has unstable states, and the largest region of instability 
in the (e^,/ 2 ) plane is obtained by choosing /3 = § and 
7 = 0. The kinetic energy of the unstable plane waves is 
bound by the condition ek < 2|<?2|h. The eigenvectors of 
the Bogoliubov matrix and the corresponding perturba- 
tions 8ij) were obtained numerically. There exists always 
two density modes Sip 1,2 , which can approximately be 
written as dip 1 ' 2 ~ Cip, where C is a time- and posi- 
tion -dependent function. The density modes are stable. 
The remaining four modes Sip 3 ' 4 ' 5,6 are approximately or- 
thogonal to ip, but it is not as easy to characterize these 
modes as in the case where the spin and magnetic field 
are parallel (/3 = 0). In general, all these modes affect 
both magnetization and spin direction. However, when 
j3 = tt/2, these modes can be classified into magnetiza- 
tion and spin modes. The magnetization mode is of the 
form (|26|) with f z = 0. This mode changes, in addition 
to the magnetization, also the spin component in the xy 
plane. The spin mode does not change the direction of 
the spin but only its amplitude /. In Fig. [3] we plot 
the positive imaginary part Wi of the eigenvalues of these 
modes. 




FIG. 3. (Color online) The amplitude of long wavelength in- 
stabilities for rubidium in the limit q ^> eu, \g2\n. Sodium 
condensate does not have long wavelength instabilities in this 
limit. The units of and u)\ are \g2\n and \g2\n/h, respec- 
tively. Now f _L B and the green color [larger lobe in (a) and 
rightmost lobe in (b)] indicates magnetization modes, while 
the blue color indicates spin modes. In (a) 7 = and in (b) 
7 = it/2. At / = 1 the figures are identical. 



V. SPIN AND MAGNETIC FIELD 
ORTHOGONAL 

In this section we compare the stability properties of 
states with f || B and f 1 B. We argue that the energies 
of unstable plane waves for states with f B are almost 
always smaller than the corresponding energies of the 
f || B case. This claim is based on energetic arguments. 
The kinetic energy ek of the domain structure can be as- 
sumed to increase as the energy of the initial state (with 
fixed magnetization) increases. The energy of the Zee- 
man term, q(l — po), is maximized when po = 0, which is 
the case if and only if the initial state is ip\\ . Furthermore, 
for a rubidium condensate also the interaction energy is 
maximized by ip» because then ginf 2 /2 = — \gi\nf 2 /2, 
which is the largest possible spin interaction energy for 
a homogeneous state with magnetization f z . For sodium 
the situation is more complicated. For q 3> gin the 
magnetic-field energy dominates and maximizes the 
energy. On the other hand, if q < gin, the energy is 
maximized when / « 1. As in the case B || f, states 
corresponding to the largest possible kinetic energy of 
the domain structure can be expected to be those with 
f z = 0. Therefore in the following we assume that mag- 
netization vanishes. It is easy to show that under this 
condition ^11 (with f z = 0) is the state with highest en- 
ergy if q > 2gin. On the other hand, when q — 0, the 
energy is maximized by 



(36) 



for which / = 1 and which is unique up to a global phase 
and a rotation around the z axis. We now compare the 
stability of this state to that of ip\\ ■ Numerically, it can be 

shown that for this state the operator Hb is periodic and 
it is therefore possible to use Floquet analysis to study 
the stability. The Floquet theorem (see, e.g., [26|) states 
that if Hb is periodic, the time evolution operator Ub 
determined by equation (fT12j) can be written as 




U B (t)=M(t)e 



itK 



(37) 



where M is a periodic matrix with period T and M(0) =1 
and K is some time-independent matrix. At times t = 
nT, where n is an integer, we get C/s(nT) = e~ mTK . The 
eigenvalues of K determine the stability of the system. 
If Ub (T) were unitary, all the eigenvalues of K would be 
real. In our case Ub{T) does not have to be unitary and 
the eigenvalues of K can have a nonvanishing imaginary 
part. We say that the system is unstable if at least one of 
the eigenvalues of K has a positive imaginary part. Wc 
denote the imaginary part of an eigenvalue uj of K by uji 
and calculate it from 



hn[i In A] 
f : 



(38) 
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where A is an eigenvalue of Ub(T). We calculated the 
eigenvalues and eigenvectors numerically for various val- 
ues of q. The oscillation period T can be obtained from 
the equations given in [271 ]. The unstable perturbations 
corresponding to the eigenvectors of K are similar to the 
ones obtained in the previous section in the (3 = n/2 
case. Hence the magnetization mode changes both mag- 
netization and the direction and length of the spin vector 
f and the spin mode affects only the length of the spin 
vector. In Fig. |4] we plot the unstable modes for some 
values of q. For comparison, also the unstable modes of 
the ip\\ states are shown. We find that for rubidium the 



(a), 



(b), 



(c), 



(d) 1 2 3 4 5 (e) 1 2 3 4 5 (f) 1 2 3 4 5 
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FIG. 4. The unstable modes of (a)-(c) rubidium and (d)- 
(f) sodium for ip± (/ = 1) and ip^ (/ = 0). Here in (a) 
and (d) q = 0.8, in (b) and (e) q = 1.2, and in (c) and (f) 
q = 1.5 in units of \g2\n. The units of £fc and uui are \g2\n and 
\g2\n/h, respectively. The dashed (dot-dashed) line gives the 
spin (magnetization) mode of i/jy , while the solid (dotted) line 
indicates the spin (magnetization) mode of ip±. 

maximal kinetic energy of the unstable perturbations of 
ip\\ is always higher than that of i(j±. For sodium the 
same conclusion holds when q > X-hg^n. If q < 1.5 (72" 
the maximal value of of can be slightly larger for ip± , 
as can be seen from Figs. HJd)-[3Jf ) . On the other hand, 
the growth rate of these instabilities is much smaller than 
the growth rate of the instabilities of ■ We therefore 
conclude that a lower bound for the wavelengths of un- 
stable perturbations is essentially given by the equation 
( k = \92\n— 92n-\-q, which is the corresponding bound for 
the states of the form tpu. Consequently, we conjecture 
that for condensate sizes smaller than the wavelength 
corresponding to = \g 2 \n — g 2 n + q both rubidium and 
sodium condensates are essentially stable regardless of 
the initial state. This wavelength is determined by 



A 



2irh 



^/2m(\g 2 \n - g 2 n + q)' 



(39) 



and wavelengths smaller than this are shown by the 
shaded region in Fig. [2] One should note that Eq. (|39|) 
gives only a sufficient condition for stability, it does not 
allow us to conclude that a condensate is unstable if it is 
larger than this size. Depending on the initial state, the 
condensate may be stable even if it is larger than the size 
determined by (|39l) . 



In addition to giving a bound for stable condensate 
size, this result makes it possible to derive constraints 
for the validity of the single- mode approximation (SMA). 
The SMA states that spatial degrees of freedom decouple 
from spin dynamics when the condensate is smaller than 
the spin healing length 



\J^rn\g 2 \ 



(40) 



This condition is obtained by requiring that the spin- 
interaction energy is insufficient to create spatial spin 
structures and its validity has been confirmed experimen- 
tally: For a 23 Na condensate with Thomas-Fermi radius 
smaller than £ s , the SMA was found to provide a very 
good description of the system [l3|. However, the va- 
lidity of SMA is also constrained by the results of the 
stability analysis discussed in this paper. If we assume 
that SMA holds initially, then the stability analysis shows 
that an additional requirement for the validity of SMA 
is that the condensate is smaller than the wavelength 
given by Eq. (|39l) . In particular, at a high magnetic field 
(<Z 3* \92\n) condition (|39|) gives a stricter bound for the 
condensate size than Eq. (|40[) . We remark that an equa- 
tion resembling Eq. (|3"9")l can be obtained also by equat- 
ing the maximal energy in a magnetic field, g 2 n/2 + q, 
with the kinetic energy . This is an extension of the ar- 
gumentation used in obtaining Eq. (j40|) to the case where 
magnetic field is nonzero. The difference between these 
approaches is that Eq (|39[) is obtained from rigorous sta- 
bility analysis, while the healing length argumentation is 
an order of magnitude estimate for the energy scales of 
the dynamics. 



VI. BOSONS ON A RING 

As a specific realization of the instabilities discussed in 
this paper we study a gas of bosonic atoms in a toroidal 
trap. We consider a doughnut-shaped condensate with 
N atoms, thickness 2p± (2p z ) in the radial (axial) direc- 
tion, and mean radius R, and assume that the trap is 
well approximated by a harmonic oscillator potential in 
the radial and axial directions, with trapping frequencies 
cjj_ and lo Z) respectively. Provided that both the spin 
healing length £ s and the wavelength given by Eq. (|3"9")l 
are larger than p± and p z , the SMA applies in radial and 
axial directions. This makes it possible to integrate out 
the dynamics in these directions. If, in addition, R is 
large enough compared to p± and p z , the condensate can 
be described as a homogeneous one-dimensional system 
of length 2irR with periodic boundaries. As a specific ex- 
ample, we discuss an optical trap of the type used in Rcf. 
[28| . created as a combination of a Laguerre-Gaussian 
beam and a laser sheet. The effective interaction energy 



32 "off 



Nh 2 (a 2 -a ) 8 
3mRp±p z 3ir ' 



(41) 
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where n e s comes from integrating the squared density 
in the Thomas-Fermi approximation in radial and axial 
directions. The Thomas-Fermi approximation can be as- 
sumed to be valid if Hu z , hioj_ -C go-off- 

As we have seen, the parameters determining insta- 
bilities are /, the angles (3 and 7, and the mode en- 
ergy in units of the spin -interaction energy, £k/\g2\n. 
In the periodic geometry considered here, k is quantized 
as k = k/R, where k is an integer. The corresponding 
mode energy is e K = H 2 k 2 /2mR 2 , and the allowed values 
for the ratio of the mode energy to the interaction energy 
are 

e* 9tt p ± p z 2 _ 2 

-. — : = : tK = e\K , 42) 

\g 2 \n cS 16 NR\a 2 -a \ ' V ' 

where for convenience we introduced the dimensionless 
prefactor ei; note that e\ = ei/|<72|^eff- The charac- 
teristic time scale for the instabilities is seen from Eqs. 
([2"T ]) -([2"5 ]) and Fig. [T] to be given by h/\g 2 \n cS (note that 
the maximum magnitude of the spin and magnetization 
modes is independent of the magnetic-field parameter 
q). With the chosen parameters, the time scale is about 
130 ms for Rb and 10 ms for Na. We simulate the time 
development of the system starting from initial states of 
the form , which we argued to be the most unstable 
ones for given magnetization. We discuss first the time 
evolution of a rubidium condensate. 
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FIG. 5. (Color online) Time development of (top) the lo- 
cal spin projection f z , (middle) squared spin / 2 , and (bot- 
tom) zeroth component po in a one-dimensional (ID) ru- 
bidium condensate with periodic boundary conditions. The 
system consists of N — 10° Rb atoms in magnetic field 
B = 130 mG (q = |<72|n e ff), initially in state ip\\ [Eq. ([17)) ] 
with / = f z — 0.2. The spatial dimension is along the verti- 
cal, and time is along the horizontal direction. 



A. Rubidium 

Figure [5] displays the time development for 87 Rb atoms 
in the initial state ip\\ with / = 0.2 and in a magnetic 
field B = 130 mG, corresponding to q = I.92 l^cff, as in 
Fig-Hfb). Note that since po and f 2 depend on squared 
wave functions, the plots exhibit second harmonics, i.e., 
the number of peaks is twice the wave number k of the 
excitation. We see that the local spin amplitude / and 
the population of the zero component po develop insta- 
bilities with dominant wave number k = 4, correspond- 
ing to £4 / 1 (72 1 ''T-cff ~ 1.91. This is close to the value 
^k/\g2\n e s — 2 which gives the fastest-growing spin mode; 
see Fig. [lib). Hence this mode is a spin mode. Another 
mode with wave number k = 3 affects both / and the 
local magnetization f z but not the zeroth spin compo- 
nent, indicating that this is a magnetization mode. For 
this mode £4 / 1 (72 1 ^eff ~ 1-07. As can be seen from Fig. 
[1Kb), this is close to the fastest-growing magnetization 
mode. We see that the linear analysis explains well the 
initial growth of the instabilities. At longer times, non- 
linear processes take over. These will not be discussed in 
more detail here. The particle density, not plotted in Fig. 
[5J stays constant to within a few percent; the instability 
only affects the spin. The time scale for buildup of an 
appreciable spin magnitude is slightly above 1 s, within 
which the modes have increased by about four orders of 
magnitude. This time is within attainable condensate 
lifetimes. 
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FIG. 6. (Color online) Time development of a ID rubidium 
condensate with periodic boundary conditions, as in Fig. [5] 
Here the initial state has a magnetization f z = f = 0.8. 



If the initial state has a higher value of /, the wave 
number and the amplitude of the most unstable magne- 
tization mode are decreased; see Fig. [T](b). An example 
for / = f z = 0.8 is given in Fig. |6] The spin mode 
still has wavenumber k = 4, which is consistent with the 
fact that the location of the fastest-growing spin mode 
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FIG. 7. (Color online) Time development of a ID rubidium 
condensate with periodic boundary conditions, as in Fig. [S] 
Here the system is made smaller so that the parameter ei = 2, 
and the external magnetic field is B = 10 mG, corresponding 
to q = |(?2|wcff. The initial state has a magnetization f z = 
/ = 0.4. 



does not depend on /. The wave number of the most 
unstable magnetization mode is reduced to k = 2. This 
gives £2/|32|«off ~ 0.48, while the fastest-growing mag- 
netization mode can be calculated from Eq. (|22|) to be 
at ek/\g2\n e B ~ 0.36. Now it takes about 4 s for the 
instability to build up. 

Assume next that the magnetic field vanishes and the 
trap parameters arc tuned so that e\ = 2. From Eq. 
(|4"2")l we see that this can be done by, e.g., loosening the 
ring trap and decreasing the number of particles. Then 
the lowest modes, located at k = and k = 1, give 
eo/l^l^off = and ei / 1 <72 l^ioff = 2. Comparison with 
Fig. 2] shows that now all states are stable. This is 
also what we see in the simulations (not shown here). 
However, by increasing the magnetic field we may once 
again make the system unstable. In Fig. [7] we report on a 
simulation where e\ = 2 and q = Ig^l^eff, corresponding 
to B = lOmG if the radius R is left unchanged. In the 
initial state f z = .f = 0.4. In such a magnetic field, we 
expect the spin mode to be the only unstable mode, with 
wave number k = 1. The increase in po caused by the 
spin mode is clearly visible in Fig. [7J An oscillation with 
wave number k = 2 is seen to develop in the magneti- 
zation simultaneously; this is not predicted by the linear 
analysis since k = 2 lies outside the unstable region in 
this case. However, a closer look at the Fourier transform 
of the spin components shows that this is not due to a 
linear instability but is a nonlinear effect. In Fig. [5J we 
see an exponential rise of the population po.i, i-e., the 
K = 1 plane wave component of the uif = spin compo- 
nent. Populations in the mp = ±1 components, both in 
K = 1 and k = 2, are excited as secondary instabilities. 



FIG. 8. Populations of the lowest plane-wave components of 
the system in Fig. [7] (top) Populations pi. K of spin component 
nriF = 1, (middle) spin component rriF = 0, and (bottom) spin 
component ttif = — 1. Solid lines show k = 1, dashed lines 
show ft = 2, and dotted lines show k = 3. 




FIG. 9. (Color online) Time development of a ID sodium con- 
densate with periodic boundary conditions, as in Fig. [5] Here 
we simulate 23 Na atoms in an initial state with magnetization 
/ = 0.4, and the magnetic field is B — 95 mG (g = g2n e g/2). 



B. Sodium 

We now consider a system of 23 Na atoms in a toroidal 
trap with the same Thomas-Fermi length parameters as 
above. For these parameters, e\ = 0.033. The system is 
stable in zero field, as seen in Fig. 1(d). If, on the other 
hand, q = g2n e s/2 (B = 95 mG) [cf. Fig. 1(e), where 
Q = .92^off]j the main instability develops at k = 0. The 
result of the simulation is shown in Fig. [SJ For this simu- 
lation we chose an initial state with / = 0.4, which allows 
instabilities with k — 0,1,2,3. Indeed, instabilities now 
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FIG. 10. Populations of the lowest plane-wave components of 
the system in Fig. [5] Panels are as in Fig. [5] Solid lines show 
K = 0, dashed lines show k = 2, dotted lines show k = 3, and 
dash-dotted lines show k = 4. 



develop with wave numbers from k = up to 3. This is 
more clearly seen in the plot of the Fourier components in 
Fig. 1101 (where the k = 1 component, whose time depen- 
dence is similar to that of the k = 2 and k = 3 compo- 
nents, is left out in order not to clutter the figure). It is 
seen that the most unstable mode has wavenumber n — 
and corresponds to uniformly populating the m = spin 
component. 

The results reported in this section indicate that the 
dynamical instabilities studied in Sees. IIVI and fVl can be 
readily studied in existing traps and that the wave num- 
ber of the unstable modes can be controlled by managing 
the system size and magnetic field. Systems small enough 
to be stable seem to be within reach. Time scales are also 
clearly tunable. 



condensate has long-wavelength instabilities which are 
independent of the strength of the magnetic field. These 
do not exist in a 23 Na condensate. Additionally, insta- 
bilities whose wavelengths depend on the strength of the 
magnetic field are possible in both systems. For rubid- 
ium these exist already at zero field, while for sodium 
nonzero magnetic field is required for the instability to 
appear. 

The stability of long wavelength perturbations was 
solved analytically also in the case where the magnetic 
-field energy is much larger than the spin interaction en- 
ergy and the kinetic energy of the plane wave perturba- 
tions. The wavelengths of the unstable long wavelength 
perturbations are bounded by the condition ek < 2\g2\n 
regardless of the initial state. 

It was also argued that states with spin parallel to the 
magnetic field are the ones whose instabilities have the 
highest energy. This claim was based on energetic argu- 
ments and a numerical study of the stability of a state 
that is orthogonal to the magnetic field. The results al- 
low us to derive an analytical formula giving a sufficient 
condition for the size of a stable condensate at a given 
magnetic field. Condensates smaller than the size given 
by Eq. (f3U| are guaranteed to be stable. However, all 
condensates larger than this are not necessarily unsta- 
ble; if prepared in a suitable state, the system may be 
stable even if it is larger than this size. Equation (f3T))) 
gives also a criterium for the validity of the single-mode 
approximation. At a high magnetic field this condition 
gives a stricter bound for the condensate size than the 
standard condition, given by the spin healing length. 

Finally, the stability properties predicted by the linear 
Bogoliubov theory were studied by solving the Gross- 
Pitaevskii equations numerically in a ID circular geom- 
etry. It was shown that by controlling the number of 
particles, trapping frequencies, and strength of the mag- 
netic field it is possible to control the stability properties 
of the condensate. 



VII. CONCLUSIONS 

We have studied the stability of spin-1 Bosc-Einstein 
condensates, concentrating on the nonstationary states 
of rubidium and sodium condensates. The analysis was 
performed in a frame of reference where the state un- 
der investigation is stationary. The stability analysis was 
done using the Bogoliubov approach, that is, expanding 
the time evolution equations of the system to first order 
with respect to a small perturbation in the stationary 
state wave function. The resulting time evolution equa- 
tions for the perturbations were solved analytically and 
numerically, assuming that the unperturbed system is 
spatially homogeneous. In particular, the effect of an ex- 
ternal homogeneous magnetic field was examined. We 
found that the eigenmodes and eigenvectors of the per- 
turbations can be determined analytically if the spin and 
magnetic field are parallel, regardless of the strength of 
the magnetic field. These eigenmodes show that a 87 Rb 
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Appendix A 

Here we examine the time evolution of spin states by 
looking at the time evolution equations of the system. 
An arbitrary spin state can be written as 



(Al) 
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Writing ip in this way and using Eq. ([5]) give the time 
evolution equations 



Here in the first equation j 
one j = 5, 6, and 



1, 2, 3, 4 and in the second 



.dpi dp 



= h- 



dt 



hdpo 
2 dt 



92np 



dt 
d6 ±1 

ft -»T = ~92n ( V(l - Po) 2 - /| cos 9 + 1 - 



V , (l-po) 2 -L 2 sin9, 



.92" ( Po J | _ ^° ^ j cos 9 + po ± fz j -q±P, 



dt 
d 

h—Q = 2g 2 n 
at 



Po) , 



_ /(.9o + .92) + Sj V (.9o - ^2) 2 + 45o32 P (vt „, 
(.90-32) V / W f 



where we have defined s\ = s 2 = —S3 = — S4 = 1. The 
corresponding perturbations become 



(l-po)(2A,-l) + / 2 



v/(1-Po) 2 -/ 2 2 

e = 200-0! -6-1. 



cos e + 2p - 1 + 2q, 



In deriving these equations we have neglected the term 
proportional to the identity operator as it changes only 
the global phase. Clearly, if po = in the initial state, the 
populations will remain constant during the subsequent 
time evolution. This means that only the phases of the 
state ^|| , given in Eq. (|T7|) . can evolve in time. Another 
special case is obtained when f z = 0, which corresponds 
to the spin vector lying in the xy plane. In this case 
0i(i) = 0_i (t) [assuming that ^(0) = 0_i(O)]. Because 
the time evolution of and po is periodic (modulo 2w) 
with the same period, also the time evolution of the state 
vector is periodic, up to a global phase. This can be seen 
by redefining the phases as 0' k (t) — 6k(t) — 0i(t), which 
gives 0[(t) = 6'_x{t) = Q,6' (t) = B(t)/2. Although the 
state vector is periodic in time, numerical calculations 
show that in general the Bogoliubov matrix Hg is not 
periodic. An exception is given by the state ipx- F° r 
this state Hb is periodic and the stability of %p± can be 
analyzed using Floquet theory. 

It is possible to obtain an approximate propagator for 
state (|Alj) under the assumption that q ^> \g2\n. Then 
Q(t) w 0(0) + 2qt/H, which leads to rapidly oscillating 
sin O and cos and we can average over one oscillation 
period, obtaining sin© m cos0 rj 0. This gives p k = 0, 
and we get the propagator 




8V=YC j F\ 



(B4) 



where 



F = 



hjjj cos(k ■ r + uijt) + ie k sin(k ■ r + Wjt), ujj E 



+ie fe )eTl^l t sin(k-r), 



Uj = ±i\uj\, 
(B5) 



and Cj is an arbitrary nonzero complex number and j 
1,2,3,4. For j = 5,6 we get 



qt/h 



g 2 n 



- (-cfc - gin + q + fiw j )e- i( - k - r+u ^] 1 . (B6) 



In order to derive an approximate expression for Sip 3 , we 
expand ctj in Taylor series with respect to g 2 /go- We get 



U 4 



„-itg2n{l-p )/h e -it{(g2nf z -p)F z + (g 2 n(2p -l)+q)F^}/h 

(A2) 



[1 + O(g 2 /g )} 



(B7) 



Appendix B: Eigenvectors 

In the case where the magnetic field and spin are paral- 
lel the eigenvectors of Hb can be calculated analytically 
and are given, up to a normalization, by 

Xj = (<x,- (e fc + Huj), 0, e k + tiUj,aj (e k - ftu>j), 0, e k - fuvj), 

(Bl) 



For rubidium and sodium I32I/50 *C 1, which allows us to 
include only the zeroth order term in the above equation. 
This gives 



5tI> j = 



jCjF 







vo^m ^. V(1 _ s . /)/2y 



= (0,g 2 n^~pe lqt/h , 0,0,0,0) 

+ (0,0,0, 0, (-e fe - g 2 n + q + hoj J )e- ,,}t/h ,0). (B2) Here j = 1, 2, 3,4. 
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